* Nicolas Van de Sijpe, 28-May-2013
* Code to replicate JDS paper "Fungibility of health aid revisited".
* Please cite the paper if you use any of this material.


*******************************************************
* SET UP
*******************************************************
clear all
set more off, perm
local Luetalcommentpath "C:\Users\Nicolasvds\Documents\My papers\Lu et al comment"
local Luetalcommentpathdata "`Luetalcommentpath'\Data"
local Luetalcommentpathresults "`Luetalcommentpath'\Results"


****************************************************************************************
* LU ET AL SETUP
****************************************************************************************
/* Data courtesy of David Roodman (http://blogs.cgdev.org/globalhealth/2012/09/the-aid-fungibility-debate-and-medical-journal-peer-review.php).
To recreate the results that follow, use the "Single-imputed data set" in that link (which I have renamed "Lu et al single imputation.dta"). 
This is the dataset used by Lu et al in their 2010 Lancet paper. */
use "`Luetalcommentpathdata'\Lu et al single imputation.dta"
d, fulln

/* Following block of code comes from Roodman: */
ren mghegdp_who_s ghegdp_who_s
ren mghegdp_john_s ghegdp_imf_s
ren pridisgdp_up_net pridisbgdp_up_net
ren mgdppc_usd06_imf gdppc_usd_imf
ren mggegdpwb ggegdp_wb
encode country, gen(countrynumber)
tsset countrynumber year
drop if country=="Angola" | country=="Eritrea"

/* Label variables: */
label variable ghegdp_who_s "Government own health spending/GDP, WHO"
label variable ghegdp_imf_s "Government own health spending/GDP, IMF"
label variable govdisgdp_up_net "DAH to government/GDP"
label variable pridisbgdp_up_net "DAH to non-government/GDP"
label variable gdppc_usd_imf "GDP per capita"
label variable ggegdp_wb "General government exp./GDP"
label variable drdisgdp "Debt relief/GDP"
d, fulln
su

***************************************************************
* SIZE DAH-Gov VERSUS DAH-NonGov
***************************************************************
/* Check size of health aid via NGOs vs. health aid to government. */
su govdisgdp_up_net 
scalar DAHgovmean = r(mean)
su pridisbgdp_up_net
scalar DAHnongovmean = r(mean)
scalar ratiogovtonongovDAH = DAHgovmean/DAHnongovmean
di "Ratio of DAH-Gov to DAH-Nongov is "ratiogovtonongovDAH
// so DAH-Gov about 4 times as large as DAH-NonGov (as mentioned in footnote 3 of my paper).


*******************************************************
* REPLICATE FE RESULTS IN LU ET AL APPENDIX TABLE A5
*******************************************************
// WHO data.
xtreg ghegdp_who_s l.ghegdp_who_s govdisgdp_up_net pridisbgdp_up_net drdisgdp gdppc_usd_imf ggegdp_wb hiv_prevalence, fe vce(cluster country)
// IMF data:
xtreg ghegdp_imf_s l.ghegdp_imf_s govdisgdp_up_net pridisbgdp_up_net drdisgdp gdppc_usd_imf ggegdp_wb hiv_prevalence, fe vce(cluster country)


capture log close
***************************************************************************
* RECOVER BETA_ON FROM "NAIVE" ESTIMATE OF FUNGIBILITY -- FE REGRESSIONS
***************************************************************************
/* This replicates table 1 in my paper (IMF data), and does the same for the WHO data (results mentioned but not reported in paper) */
log using "`Luetalcommentpathresults'\Lancetcomment table beta.log", replace
foreach dataset in who imf { // loop over the two datasets to check whether results differ or not.
	// run FE regressions, pick up coefficient and standard error.
	xtreg ghegdp_`dataset'_s l.ghegdp_`dataset'_s govdisgdp_up_net pridisbgdp_up_net drdisgdp gdppc_usd_imf ggegdp_wb hiv_prevalence, fe vce(cluster country)
	scalar beta_`dataset' = _b[govdisgdp_up_net]
	scalar se_beta_`dataset' = _se[govdisgdp_up_net]
	// this is the "naive" estimate, plus its standard error.
	// also pick up coefficient on lagged dependent variable.
	scalar lagdepcoeff`dataset' = _b[L.ghegdp_`dataset'_s]
	// also pick up the LR effect:
	nlcom LReffect: _b[govdisgdp_up_net]/(1-_b[L.ghegdp_`dataset'_s]), post
	scalar define betaLR_`dataset' = _b[LReffect]
	scalar define se_betaLR_`dataset' = _se[LReffect]

	
	// generate matrices of the appropriate size, containing 0's.
	foreach matrix in beta_on_`dataset' beta_on_`dataset'_se beta_on_`dataset'_fullfungp beta_on_`dataset'_nofungp ///
	betaLR_on_`dataset' betaLR_on_`dataset'_se betaLR_on_`dataset'_fullfungp betaLR_on_`dataset'_nofungp {
		matrix `matrix' = J(7,7,0)
		matrix colnames `matrix' = -3/4 -1/2 -1/4 0 1/4 1/2 3/4 // columns = correlation between on- and off-budget aid (rho).
		matrix rownames `matrix' = 1/4 1/2 3/4 1 3/2 2 4 // rows = relative variance of off-budget versus on-budget aid (delta).
	} // end generation zero matrices

	qui {
		local i = 0
		foreach relvar_off_vs_on of numlist 0.25 0.5 0.75 1 1.5 2 4 { // vary relative variance of off-budget versus on-budget aid (delta).
			local++i
			local relstdev_off_vs_on = sqrt(`relvar_off_vs_on')
			local j = 0
			foreach corr_on_off of numlist -0.75 -0.5 -0.25 0 0.25 0.5 0.75 { // vary correlation between on- and off-budget aid (rho).
				local ++j
				// SHORT-RUN EFFECTS FIRST:
				/* coefficient "beta_on" entry for matrix */
				scalar entry`i'`j'_beta_on_`dataset' = beta_`dataset' + 1 - ((1 + `corr_on_off' * `relstdev_off_vs_on')/(1 + `relstdev_off_vs_on'^2 + 2 * `corr_on_off' * `relstdev_off_vs_on')) 
				/* standard error is just the standard error of the coefficient in the original regression. */
				scalar entry`i'`j'_se_`dataset' = se_beta_`dataset'
				/* Re-run FE regression: */
				xtreg ghegdp_`dataset'_s l.ghegdp_`dataset'_s govdisgdp_up_net pridisbgdp_up_net drdisgdp gdppc_usd_imf ggegdp_wb hiv_prevalence, fe vce(cluster country)
				/* full fungibility test (H_0 is coefficient smaller than or equal to -1): */
				test _b[govdisgdp_up_net] = -1
				// calculate p-value for one-sided full fungibility test:
				scalar fullfungtest`i'`j'_`dataset' = ttail(r(df_r),(entry`i'`j'_beta_on_`dataset'+1)/entry`i'`j'_se_`dataset')
				// no fungibility test (H_0 is coefficient larger than or equal to 0)
				test _b[govdisgdp_up_net] = 0
				// calculate p-value for one-sided no fungibility test:
				scalar nofungtest`i'`j'_`dataset' = 1 - ttail(r(df_r),entry`i'`j'_beta_on_`dataset'/entry`i'`j'_se_`dataset')

				// NOW LONG-RUN EFFECTS:
				/* LR effect of DAH-Gov first. */
				scalar entry`i'`j'_betaLR_on_`dataset' = entry`i'`j'_beta_on_`dataset'/(1-lagdepcoeff`dataset')
				/* standard error is just the standard error of the LR effect in the original regression. */
				scalar entry`i'`j'_seLR_`dataset' = se_betaLR_`dataset'
				/* post LR effect: */
				nlcom LReffect: _b[govdisgdp_up_net]/(1-_b[L.ghegdp_`dataset'_s]), post // this overrides previous estimation results
				/* full fungibility test (H_0 is coefficient smaller than or equal to -1): */
				test _b[LReffect] = -1
				// calculate p-value for one-sided full fungibility test:
				scalar fullfungtestLR`i'`j'_`dataset' = ttail(r(df_r),(entry`i'`j'_betaLR_on_`dataset'+1)/entry`i'`j'_seLR_`dataset')
				// no fungibility test (H_0 is coefficient larger than or equal to 0)
				test _b[LReffect] = 0
				// calculate p-value for one-sided no fungibility test:
				scalar nofungtestLR`i'`j'_`dataset' = 1 - ttail(r(df_r),entry`i'`j'_betaLR_on_`dataset'/entry`i'`j'_seLR_`dataset')
				
				/* Build up the matrices. */
				// SR effects:
				mat beta_on_`dataset'[`i',`j'] = entry`i'`j'_beta_on_`dataset' // matrix with "beta_ON" coefficients.
				mat beta_on_`dataset'_se[`i',`j'] = entry`i'`j'_se_`dataset' // matrix with standard errors.
				mat beta_on_`dataset'_fullfungp[`i',`j'] = fullfungtest`i'`j'_`dataset' // matrix with p-values full fungibility test.
				mat beta_on_`dataset'_nofungp[`i',`j'] = nofungtest`i'`j'_`dataset' // matrix with p-values no fungibility test.
				// LR effects:
				mat betaLR_on_`dataset'[`i',`j'] = entry`i'`j'_betaLR_on_`dataset' // matrix with LR "beta_ON" coefficients.
				mat betaLR_on_`dataset'_se[`i',`j'] = entry`i'`j'_seLR_`dataset' // matrix with LR standard errors.
				mat betaLR_on_`dataset'_fullfungp[`i',`j'] = fullfungtestLR`i'`j'_`dataset' // matrix with p-values LR full fungibility test.
				mat betaLR_on_`dataset'_nofungp[`i',`j'] = nofungtestLR`i'`j'_`dataset' // matrix with p-values LR no fungibility test.
			}/* end loop over corr_on_off*/
		}/* end loop over relstdev_off_vs_on */
	}/* end qui */
	di "SR EFFECTS DATASET `dataset'"
	matrix list beta_on_`dataset'
	matrix list beta_on_`dataset'_se
	matrix list beta_on_`dataset'_fullfungp
	matrix list beta_on_`dataset'_nofungp
	di "LR EFFECTS DATASET `dataset'"
	matrix list betaLR_on_`dataset'
	matrix list betaLR_on_`dataset'_se
	matrix list betaLR_on_`dataset'_fullfungp
	matrix list betaLR_on_`dataset'_nofungp
}/* end loop over datasets. */
log close
// EXPORT MATRICES
// SR effect WHO dataset:
outtable using "`Luetalcommentpath'\matrix_who_FE_SR" ///
, mat(beta_on_who) format(%4.2f) replace center nobox caption("Short-run fungibility of DAH-Gov (WHO data)") clabel("Table_matrix_varying_on_off_budget_aid_who_SR")
// LR effect WHO dataset:
outtable using "`Luetalcommentpath'\matrix_who_FE_LR" ///
, mat(betaLR_on_who) format(%4.2f) replace center nobox caption("Long-run fungibility of DAH-Gov (WHO data)") clabel("Table_matrix_varying_on_off_budget_aid_who_LR")
// SR effect IMF dataset:
outtable using "`Luetalcommentpath'\matrix_imf_FE_SR" ///
, mat(beta_on_imf) format(%4.2f) replace center nobox caption("Short-run fungibility of DAH-Gov (IMF data)") clabel("Table_matrix_varying_on_off_budget_aid_imf_SR")
// LR effect IMF dataset:
outtable using "`Luetalcommentpath'\matrix_imf_FE_LR" ///
, mat(betaLR_on_imf) format(%4.2f) replace center nobox caption("Long-run fungibility of DAH-Gov (IMF data)") clabel("Table_matrix_varying_on_off_budget_aid_imf_LR")


capture log close
********************************************************************************************************
* RECOVER BETA_ON FROM "NAIVE" ESTIMATE OF FUNGIBILITY STARTING FROM BIAS-CORRECTED LSDV REGRESSIONS
********************************************************************************************************
/* This replicates table 2 in my paper (IMF data), and does the same for the WHO data (results mentioned but not reported in paper) */
log using "`Luetalcommentpathresults'\Lancetcomment table beta LSDVC.log", replace
clear results 
clear matrix
local bootstrapreps 100 // number of boostrap replications used to construct standard errors.
foreach dataset in who imf { // loop over the two datasets just to check whether results are different.
	// run LSDVC regressions, pick up coefficient and standard error, plus store results for later.
	eststo `dataset'results: xtlsdvc ghegdp_`dataset'_s /*l.ghegdp_`dataset'_s*/ govdisgdp_up_net pridisbgdp_up_net drdisgdp gdppc_usd_imf ggegdp_wb /// 
	hiv_prevalence, initial(bb) bias(1) vcov(`bootstrapreps')
	scalar beta_`dataset' = _b[govdisgdp_up_net]
	scalar se_beta_`dataset' = _se[govdisgdp_up_net]
	// this is the "naive" estimate, plus its standard error.
	// also pick up coefficient on lagged dependent variable.
	scalar lagdepcoeff`dataset' = _b[L.ghegdp_`dataset'_s]
	// also pick up the LR effects:
	nlcom LReffect: _b[govdisgdp_up_net]/(1-_b[L.ghegdp_`dataset'_s]), post
	scalar define betaLR_`dataset' = _b[LReffect]
	scalar define se_betaLR_`dataset' = _se[LReffect]

	
	// generate matrices of the appropriate size, containing 0's.
	foreach matrix in beta_on_`dataset' beta_on_`dataset'_se beta_on_`dataset'_fullfungp beta_on_`dataset'_nofungp ///
	betaLR_on_`dataset' betaLR_on_`dataset'_se betaLR_on_`dataset'_fullfungp betaLR_on_`dataset'_nofungp {
		matrix `matrix' = J(7,7,0)
		matrix colnames `matrix' = -3/4 -1/2 -1/4 0 1/4 1/2 3/4 // columns = correlation between on- and off-budget aid.
		matrix rownames `matrix' = 1/4 1/2 3/4 1 3/2 2 4 // rows = relative variance of off-budget versus on-budget aid.
	} // end generation zero matrices

	qui {
		local i = 0
		foreach relvar_off_vs_on of numlist 0.25 0.5 0.75 1 1.5 2 4 { // vary relative variance of off-budget versus on-budget aid.
			local++i
			local relstdev_off_vs_on = sqrt(`relvar_off_vs_on')
			local j = 0
			foreach corr_on_off of numlist -0.75 -0.5 -0.25 0 0.25 0.5 0.75 { // vary correlation between on- and off-budget aid.
				local ++j
				// SHORT-RUN EFFECTS FIRST:
				/* coefficient "beta_on" entry for matrix */
				scalar entry`i'`j'_beta_on_`dataset' = beta_`dataset' + 1 - ((1 + `corr_on_off' * `relstdev_off_vs_on')/(1 + `relstdev_off_vs_on'^2 + 2 * `corr_on_off' * `relstdev_off_vs_on')) 
				/* standard error is just the standard error of the coefficient in the original regression. */
				scalar entry`i'`j'_se_`dataset' = se_beta_`dataset'
				/* Restore LSDVC regression: */
				estimates restore `dataset'results
				/* full fungibility test (H_0 is coefficient smaller than or equal to -1): */
				test _b[govdisgdp_up_net] = -1
				// calculate p-value for one-sided full fungibility test:
				scalar fullfungtest`i'`j'_`dataset' = 1 - normal((entry`i'`j'_beta_on_`dataset'+1)/entry`i'`j'_se_`dataset')
				// no fungibility test (H_0 is coefficient larger than or equal to 0)
				test _b[govdisgdp_up_net] = 0
				// calculate p-value for one-sided no fungibility test:
				scalar nofungtest`i'`j'_`dataset' = normal(entry`i'`j'_beta_on_`dataset'/entry`i'`j'_se_`dataset')

				// NOW LONG-RUN EFFECTS:
				/* LR effect of DAH-Gov first. */
				scalar entry`i'`j'_betaLR_on_`dataset' = entry`i'`j'_beta_on_`dataset'/(1-lagdepcoeff`dataset')
				/* standard error is just the standard error of the LR effect in the original regression. */
				scalar entry`i'`j'_seLR_`dataset' = se_betaLR_`dataset'
				/* post LR effect: */
				nlcom LReffect: _b[govdisgdp_up_net]/(1-_b[L.ghegdp_`dataset'_s]), post // this overrides previous estimation results
				/* full fungibility test (H_0 is coefficient smaller than or equal to -1): */
				test _b[LReffect] = -1
				// calculate p-value for one-sided full fungibility test:
				scalar fullfungtestLR`i'`j'_`dataset' = 1 - normal((entry`i'`j'_betaLR_on_`dataset'+1)/entry`i'`j'_seLR_`dataset')
				// no fungibility test (H_0 is coefficient larger than or equal to 0)
				test _b[LReffect] = 0
				// calculate p-value for one-sided no fungibility test:
				scalar nofungtestLR`i'`j'_`dataset' = normal(entry`i'`j'_betaLR_on_`dataset'/entry`i'`j'_seLR_`dataset')
				
				/* Build up the matrices. */
				// SR effects:
				mat beta_on_`dataset'[`i',`j'] = entry`i'`j'_beta_on_`dataset' // matrix with "beta_ON" coefficients.
				mat beta_on_`dataset'_se[`i',`j'] = entry`i'`j'_se_`dataset' // matrix with standard errors.
				mat beta_on_`dataset'_fullfungp[`i',`j'] = fullfungtest`i'`j'_`dataset' // matrix with p-values full fungibility test.
				mat beta_on_`dataset'_nofungp[`i',`j'] = nofungtest`i'`j'_`dataset' // matrix with p-values no fungibility test.
				// LR effects:
				mat betaLR_on_`dataset'[`i',`j'] = entry`i'`j'_betaLR_on_`dataset' // matrix with LR "beta_ON" coefficients.
				mat betaLR_on_`dataset'_se[`i',`j'] = entry`i'`j'_seLR_`dataset' // matrix with LR standard errors.
				mat betaLR_on_`dataset'_fullfungp[`i',`j'] = fullfungtestLR`i'`j'_`dataset' // matrix with p-values LR full fungibility test.
				mat betaLR_on_`dataset'_nofungp[`i',`j'] = nofungtestLR`i'`j'_`dataset' // matrix with p-values LR no fungibility test.
			}/* end loop over corr_on_off*/
		}/* end loop over relstdev_off_vs_on */
	}/* end qui */
	di "SR EFFECTS DATASET `dataset'"
	matrix list beta_on_`dataset'
	matrix list beta_on_`dataset'_se
	matrix list beta_on_`dataset'_fullfungp
	matrix list beta_on_`dataset'_nofungp
	di "LR EFFECTS DATASET `dataset'"
	matrix list betaLR_on_`dataset'
	matrix list betaLR_on_`dataset'_se
	matrix list betaLR_on_`dataset'_fullfungp
	matrix list betaLR_on_`dataset'_nofungp
}/* end loop over datasets. */
log close
// EXPORT MATRICES
// SR effect WHO dataset:
outtable using "`Luetalcommentpath'\matrix_who_LSDVC_SR" ///
, mat(beta_on_who) format(%4.2f) replace center nobox caption("Short-run fungibility of DAH-Gov (WHO data), bias-corrected estimator") /// 
clabel("Table_matrix_varying_on_off_budget_aid_who_SR_LSDVC")
// LR effect WHO dataset:
outtable using "`Luetalcommentpath'\matrix_who_LSDVC_LR" ///
, mat(betaLR_on_who) format(%4.2f) replace center nobox caption("Long-run fungibility of DAH-Gov (WHO data), bias-corrected estimator") /// 
clabel("Table_matrix_varying_on_off_budget_aid_who_LR_LSDVC")
// SR effect IMF dataset:
outtable using "`Luetalcommentpath'\matrix_imf_LSDVC_SR" ///
, mat(beta_on_imf) format(%4.2f) replace center nobox caption("Short-run fungibility of DAH-Gov (IMF data), bias-corrected estimator") /// 
clabel("Table_matrix_varying_on_off_budget_aid_imf_SR_LSDVC")
// LR effect IMF dataset:
outtable using "`Luetalcommentpath'\matrix_imf_LSDVC_LR" ///
, mat(betaLR_on_imf) format(%4.2f) replace center nobox caption("Long-run fungibility of DAH-Gov (IMF data), bias-corrected estimator") /// 
clabel("Table_matrix_varying_on_off_budget_aid_imf_LR_LSDVC")
